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ABSTRACT 


A method is developed which allows us to 
study the evolution of rotating stars well beyond 
the main sequence stage. We consider four different 
cases of redistribution of angular momentum in an 
evolving star. Evolutionary sequences for a 7M^ 
star, rotating according to these different cases, 
were computed from the ZAAJS to the double shell 
source stage. Each sequence was begun with a 
(typical) equatorial rotational velocity of 210 km/sec. 

On the main sequence, the effects of rotation are 

of minor importance. However, as the core contracts / 

during later stages , important effects arise in , / 

all physically plausible cases. The outer regions 

of the cores approach critical velocities and develop 

unstable angular velocity distributions. The effects 

of these instabilities should significantly alter 

the subsequent evolution. 
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I. INTRODUCTION 

In the past twenty years, astonishing progress 
has been mads toward understanding stellar 
evolution. In spite of an enormous number of 
studies of spherical models, however, relatively 
little attention has been paid to the role of 
rotation. Basically, there are two reasons for 
the neglect of rotation in most previous investi- 
gations: first, dropping the assumption of spherical 

symmetry leads to a substantial increase in the 
numerical complexity of the equations of stellar 
structure and, second, spherically-symraetric (non- 
rotating) models have been very successful in 
explaining the relevant observational data (the 
mass-luminosity relationship, H-R diagrams of 
clusters , etc . ) . 

The problem of numerical complexity can be 
reduced somewhat by using coordinate systems 
defined in terms of equipotential surfaces. This 
will be discussed further in Section II. The 
success of spherical stellar evolution theory can 
be understood by considering the stages of 
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evolution to which most observational tests apply, 
i.e., the main sequence (MS) and early post-WS 
stages. Calculations of rotating MS stars , for 
example, show that the effects of rotation on 
the internal structure of such stars are very 
small, unless the interiors are rotating much more 
rapidly than the surface layers (see , e.g., Sackraann 
and Anand 1970, Bodenheimor 1971). Because of 
the long duration of the MS stage, it is reasonable 
to expect that various dissipation mechanisms 
(such as large scale circulation currents) will 
reduce any rapid core rotation built up during 
contraction to the MS, so rotation will not play 
a significant role during the MS and early post-MS 
stages. 

The situation is very different for the later 
stages of evolution. Consider the ratio of the 
centrifugal force, , due to rotation to the 
force, F , due to gravity in the equatorial 

O 

plane. As successive fuels are exhausted and the 
core contracts, will increase in inverse 

proportion to the radius of the coi-e, if anguiar 
momentum is conserved within the core. Redistribution 
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of angular momentum throughout the star is inhibited 
by the short time scales for the evolution and by 
large radial gradients in the mean molecular weight, 
which choke off large scale circulation currents. 
Sackmahn and Weidemann (1972) and Maeder (1974) 
have shown that rotation should play a major role 
in the evolution subsequent to core helium exhaustion 
and Sofia (1971) has shown that the rotation rate 
may be the central parameter in the final gravitat- 
ional collapse of a star. 

Aside from calculations of differentially 
rotating white dwarfs, which we will not consider 
since they are not evolutionary sequences , we are 
aware of only two investigations in which models 
were calculated well beyond the MS with the effects 
of rotation included. Kippenhahn , Meyer-Hofmeister , 
and Thomas (1969, subsequently referred to as KMT) 

followed the evolution of a 9M star from the zero- 

o 

age main sequence (ZAMS) to the end of helium 
burning in the core. Two sequences of models were 
computed, corresponding to two prescriptions of 
how u (M^) , the angular velocity as a function of 
the mass coordinate, evolves with time. Meyer- 



Hofmeister (1972) computed sequences for 5, 6, and 
9M^ stars through the same stages of evolution but 
using only one of the prescriptions for 'v(M^) used 
by KMT. In both investigations, spherical symmetry 
was assumed and the effects of rotation were 
approximated by taking the radial component of the 
centrifugal force (averaged over a sphere) into 
account in the hydi'ostatic equilibrium equation. 

As will be discussed in Section II, this provides 
only a rough estimate of the effects of rotation 
on the evolution. 

It appears, therefore, that considerable work 
remains to be done on the role of rotation in the 
post-MS stages of evolution. This is the first 
of a series of papers in which we will systematically 
explore the effects of rotation on these stages. 

In this paper , we will concentrate on the computational 
techniques and describe the i-esults of some sample 
calculations on the evolution of a 7M^ star 
from the ZAMS to the double shell source stage, 
using several different prescriptions for the 
evolution of jli(M^). In Section II, we outline 
the method we use to include the, effects of 




rotation. Comparisons with other methods and 
numerical tests are described in Section III. 

The starting models and the prescriptions for 
tw(M^) are specified in Section IV. The results of 
the evolutionary calculations are presented in 
Section V and discussed in Section VI. A detailed 
description of the evaluation of the total potential 
is given in appendices. 

II. THE EQUATIONS FOR ROTATING STARS 
There are four ways in which rotation may 
effect the equations of stellar structure : 

(1) centrifugal forces reduce the effective 
gravity at any point , not on the axis of 
rotation. This must be taken directly into 
account in the equation of hydrostatic equil- 
ibrium. ~ 

(2) because the centrifugal foi'ce is not, in 
general, parallel to the force of gravity, 
equipotential surfaces are no longer spheres 
and the spherical relationships between the 
radius, enclosed volume, and surface area of 
an equipotential surface cannot be used. 
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This affects all of the equations except the 
equation for adiabatic convective equilibrium. 

(3) because the radiative flux varies with the 
local effective gravity (the von Zeipel effect), 
the radiative flux is not constant on an 
equipotential surface. This enters directly 
into the radiative equilibrium equation and 

may affect the stability to convection by 
changing the radia+ive temperature gradient. 

(4) rotation may inhibit cex'tain modes of convective 
motions and, thus, directly affect the criterion 
for convective stability (cf. Randei'S 194 2, 
Cowling 195J ) . In addition, not all angular 
momentum distributions are stable and this 

can lead to convection in regions which are 
stable to purely thermal convection (cf. 
Wasiutynski 1946). 

The first three effects can be incorporated 
into the equations of stellar structure in a 
fairly direct manner, if the total (gravitational 
plus rotational) potential, ii , is conservative. 

The formulation we use has been derived by Kippenhahn 
and Thomas (1970) and we will subsequently refer 
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to this formulation as the KT method. Because 
the form in which we use the equations is slightly 
different from that of Kippenhahn and Thomas and 
in order to clarify how the effects of rotation 
are taken into account , we rederive the equations 
here . 


The spherical surfaces normally used in stellar 
models are replaced by equipotential surfaces. 

The area of such a surface is denoted by S. and 


the volume enclosed by the surface by V 


V 


On 


such a surface, the pressure, P, the temperature, 

T, and the density, p, are constant, if the total 
potential is conservative. The Lagrangian coordinate, 
M , is replaced by M. , the mass interior to the 

r , if 

equipotential surface and the spatial variable, r, 
is replaced by r , the radius of a sphere enclosing 


a volume V, , i.e . , 


w - 4” 3 


( 1 ) 


For any quantity f, which is not constant over an 
equipotential surface, a mean value is defined 
by 


" S“ J ^ ^ ’ ’ 

'' i); = const 


( 2 ) 


page JS 
^OOR QUALITY 



where da is an element of the surface tji constant. 
The local effective gravity is defined by 

g = diji /dn (3) 

where dn is the distance between the surfaces 
\i| “ constant and r d,;) = constant. Because of 
the definition of r, , the form of the mass continuity 
equation is not altered by rotation: 

or 

5r^,/^)M^ = l/4nr^_^p. (5) 

From equation (2), 

c'V = j dnd-7 = a i 

’ y “ Const. “ V = const. 

(dn/d\!i)dj = g~^ S , d'li (6) 

and combining this with equation (4) gives 

= (dV^/aor^ SMVo 

aM /g~^ S 0 . (7) 

. . . . . 

V 

Equation (7) can be combined with the general form 
of the equation of hydrostatic equilibrium, 

3P/dS- = -p , (8) 


to give 
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gives , with equation (2) , 


L = 


4acT^ 

“3T~ 



( 14 ) 


Using equation (9) , we rewrite equation (14) as 


In T _ 3 h P 

a' 'in p ~ rr ¥" I" ' 

T P 

where 



V 


2 2 


1 



g g 


(15) 


(16) 


In the present formulation, we neglect the fourth 
effect of rotation on the equations (s?e above) 
and use the Schwarzschild criterion for convection- 
Then 


5 In T 

T-nrr 


min . 



^rad 



(17) 


whei-e V , and v are the normal (spherical) 
adiabatic and radiative gradients. We note that 
equation (11) is an approximation in that expansion 
and contraction will not, in general, maintain a 
conservative potential (cf. Kippenhahn and 
Mollenhoff 1974). As a result, the last two terms 
in equation (11) will not be constant over an 
equipotential surface. We assume that this effect 
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will be compensated by large scale circulation 
curx’ents. Aside fi'om this approximation and the 
neglect of the direct' effect of rotation on 
convection, no other approximations have been 
made in transforming the spherical equations to 
equations with rotation. 

Before describing the evaluation of the 
potential, a few comments on the above equations 
are in order. The first three effects of rotation, 
described at the beginning of this section, are 
contained in: the factor fp, the interpretation 

of r, thi’ough equation (1) and the evaluation of 
S, and V, , and the factor f_,, respectively. As 

ij/ l[( T 

the rate of rotation goes to zero, f 1 , f _ - 1 , 

p 1 

and r, - r so the equations reduce to their spherical 
counterparts. The ratio f,j,/fp enters into the 
radiative equilibrium equation because we evaluate 

a In T /3 In P = In T/:^M, )/(:^ In P/^M, ) , rather 

■? 

than In T/.^M, . The approximation used by K.OT 

V . ■ ■ 

and Meyer-Hofmeister (1972) is roughly equivalent 
to considering only the factor fp in equation (9). 

In general, fp deviates from unity by a much larger 
amount than does f ,j, , so this would appear to be a 
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consistent level of approximation. However, the 

ratio f„/f ~ 1/f , which must appear in the 
T p p 

criterion for convection [see equation (17)], 

deviates from unity by an amount comparable 

to the deviation of f and this should be taken 

P ■ 

into account at the level of approximation used 
by KMT and Meyer-Hofmeister. 

A detailed description of the evaluation of 
the potential and the factors f and f_, is given 

^ X 

in the Appendices . Basically, the potential is 
divided into three parts; 'jig. the spherically 
sy.mmetric part of the gravitational potential , 

the cylindrically symmetric potential directly 
due to rotation, and the cylindrically symmetric 

part of the gravitational potential due to distortion 
of the figure of the star. The evaluation of the 

X Xx 5 C C VVO I^UXCS x5 CXXVXUX. x no WAXXXM puxi. X£> 

evaluated, at a given point, as if that point 
were on the surface of a polytrope rotating with 
constant li’ and with the same ratio of mean 
interior density to centi’al density as displayed 
by the model at the point in question. The 
distortion of the mass exterior to the equipotential 




surface which contains the point is ignored. As 
we will show in Section III, this provides a very 
good approximation to the actual potential. 

Equations (5) , (9) , (11) and (17) are 
similar enough to the spherical equations that they 
can easily be incorporated into existing stellar 
evolution codes. In our case, we have used the 
Paczynski code with modifications to the input 
physics as described by Endal (1975). 

III. COMPARISON TO OTHER METHODS AND 
NUMERICAL TESTS 

Several different methods have been used in 
the past to investigate the effects of rotation. 
Before describing the numerical tests of our 
program, we will discuss the similarities and 
differences among the various methods. We will 
not include perturbation techniques such as those 
of Chandrasekhar (1933) and Sweet and Roy (1953) 
which have only been applied to polytropes or 
other highly simplified stellar models. 

D ouble -approximation method - This method has been 
used by a number of investigators (Roxburgh, 
Griffith, and Sweet 1965; Faulkner, Roxburgh, and 
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Strittmatter 1968; Strittmatter , Robertson, and 
Faulkner 1970; Sackmann and Anand 1970; Sackmann 
1970) , each of whom have introduced minor modifica- 
tions. The method consists of dividing the star 
into two parts : a core which is assumed to be 

rotating slowly in comparison to the local critical 
velocity and an envelope which contains a negligible 
amount of mass. First-order expansions in a 
rotation parameter are used in the core and the 
Laplace equation is used for the potential in the 
envelope (the mass of the envelope does not 
contribute to tho potential). For uniform rotation, 
this is a valid approximation but for differential 
rotation (rapidly rotating cores) the core 

approximation breaks down. 

2 2 
J method - a detailed description of the J 

method has been given by Papoloizou and Whelan 

(1973). It has also been used by Whelan, Papoloizou 

and Smith (1971 ), Whelan (1972), and htoss (1973). 

2 

The stellar structure equai 'cns used in the J 
method are formally equivalent to those used in 
the KT method so the only difference is in the 

2 ' 

evaluation of the total potential, . In the J 



method, the Roche approximation has generally been 
used. This is equivalent to neglecting the 
contribution of the quantity discussed in 
Section II. The Roche approximation should break 
down in regions where the ratio of mean irV;’rior 
density to central density (p/p^) and the ratio 
of angular velocity to critical angulaj velocity 
(uj/iicj,) are both significant compared to unity. 

This will be the case for rapid differential 
rotation in the core. 

SCF method - the self -consistent -field (SCF) 

method , as described by ©striker and Mark (1968), 

is primarily a method for obtaining accurate 

solutions of the total potential and hydrostatic 

equilibrium equations. The SCF method was combined 

with the complete stellar structure equations by 

Jackson (1970, see also Mark 1968). Jackson's 

method was applied to differentially rotating stars 

on the upper main sequence by Bodenheimer (1971). 

In this method the stellar structure equations 

are evaluated on equipotential surfaces, as in the 
2 

KT and J methods. The primary difference is that 
the potential equation is solved much more 
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accurately and the method is, as a result, very 
time-consuming. This makes it unsuitable for 
evolutionary calculations but it remains a very 
powerful method for testing more approximate 
techniques. 

Papaloizou and Whelan (1973) have made extensive 

comparisons of the results obtained with the various 

methods (including the KT method) for uniformly 

rotating ZAMS models. Except at low masses 

(M ^ they find that all the methods described 

above produce essentially the same results. At 

low masses the comparisons seem to be complicated 

by differences in the chemical composition and 

convective mixing lengths used in the models. 

Since we will not be computing models of low 

mass stars, we will not worry about this problem. 

Above IM^ , the reduction in the luminosity produced 

by uniform rotation is confined to <10% and this 

2 

I'eduction varies linearly with ;ju . This indicates 
that uniform rotation may be considered as a small 
perturbation on non-rotating models and, therefore, 
does not constitute a very stringent test for 
methods which will be used when rotation introduces 
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large changes in the models. Unfortunately, there 
are very few calculations with which we can 
compare our results for differentially rotating 
models . 

Bodenheimer (1971) has calculated ZAMS models 
of 15, 30, and 60M^ stars with strong differential 
rotation using the SCF method . In 
the models computed by Bodenheimer, the velocity 
was constant on cylinders , whereas we have 
specified constant angular velocities on equi- 
potential surfaces, so the results obtained by 
the two methods may not be directly comparable. 
However, Bodenheimer found that, to a good approx- 
imation, the effects of rotation on the central 
temperature (T^) , central density (pj,)i and 
luminosity (L) do not depend on the distribution 
of angular momentum but only on the total angular 
momentum (J) of the models. If we confine ourselve 
to these parameters, the comparison may be valid. 
Figure 1 shows the variations of these parameters, 
for a 30^^ star, as a function of log J. The 
SCF results are indicated by the symbols used 
in Figure 5 of Bodenheimer. The different symbols 



refer to different internal distributions of the 

angular momentum. The results obtained with our 

method are indicated by dashed lines. We have 

specified that (at the equator) be constant 

throughout our models. The values of w/i 

cr 

as various values of log J are indicated at the 
top of the figure. Up to “ 0.8, the 

differences between our models and the SCF models 
with the same value of J are generally comparable 
to the scatter introduced by different angular 
momentum distributions in the SCF models. Beyond 
'o)/uj = 0.8, we were not able to produce converged 

Cl 

models because of the large differences between 

the (non-rotating) Schwarzschild model used as an 

initial guess and the rotating models to be 

calculated by relaxation. We should emphasize 

here that the convergence problem is due to a 

poor initial guess, rather than to limitations of 

our method. In the evolutionary sequences presented 

in Section V, u)/m was as large as 0.99. In 

cr 

these calculations, there were no problems with 
convergence because the initial guesses were based 
on extrapolations of models which were already 
rapidly rotating. 
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IV. THE EVOLUTION OF A 7M^ STAR: STARTING 

MODELS AND ROTATION LAWS 

In order to explore the effects of rotation 
on the post-MS stages of evolution, five sequences 
of models were computed for a 7M^ star : one 

sequence for a non-rotating star (Case O) and 
four sequences for rotating stars with different 
assumptions about the redistribution of angular 
momentum in an evolving star (Cases 1-4). For 
each case , the sequence was begun on the ZAMS 
with a chemical composition of X = 0.7, Z = 0.03. 
The prescriptions for u(M^) used in the various 
cases are summarized in Table 1. 

In all cases, the angular momentum was assumed 
to be constant on equipotentia L surf aces , rather 
than on cylinders , as in the models of Bodenheiraer 
(1971). This has already been discussed in Section 
III and will be discussed further in Section VI. 
Also, in conserving angular momentum, deformation 
of the equipotentia 1 surfaces was ignored, i.e., 
the moment of inertia of a given mass shell was 
assumed to be that of a thin spherical shell of 
radius r^ . Cases 1-3 were begun with solid body 
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rotation and u =8.8 x 10“ Vsec . This corresponds 

to the average x'otation rate for spectral type B5 

main sequence stars , according to Abt and Hunter 

(1962). Unlike the models computed by KMT and 

Meyer-Hofmeister , our models are hot rapid rotators 

on the MS. The Case 4 sequence was begun with 

-5 

solid body rotation at ;u = 8.8 x 10 /sec in the 

regions outside of the convective core. Within 

-5 

the convective core, we set = min. [8.8 .x 10 /sec 
2 

(r /r) , 0.9 ij 1 , where r is the radius of the 
cc cr cc 

convective core. Limiting u) to 90% of was 

necessary to avoid the supercritical velocities 

•2 

near the center implied by the or = constant 
prescription. This limit on o was also applied 
to the convective regions of the evolving models 
for Case 4. The extra angular momentum was 
uniformly redistributed throughout the remainder 
of the convective region whenever the limiting 
velocity was reached. 

The rotation laws chosen are not meant to 
accurately reflect the effects of the various 
angulax' momentum redistribution mechanisms 
which may operate in an evolving star. Rather, 
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they are meant to bracket the set of physically 
plausible rotation laws. Cases 1 and 3 represent 
extremes in terms of complete redistribution of 
angular momentum and no redistribution of angular 
momentum, respectively. Case 2 represents an 
intermediate (and more likely) possibility that, 
while circulation currents redistribute the 
angular momentum in chemically homogeneous radiative 
regions, such currents are choked off in inhomogeneous 
regions by gradients in the mean molecular weight 
(cf. Mestel 1953). In convective regions, the 
high viscosity associated with turbulence may 
lead to solid body rotation. Case A represents 
an alternate possibility, namely that mass motions 
in convective regions tend to equalize the specific 
angular momentum in such regions. This would be 
the case if the circulation produced by convection . , 
preserves the angular momentum of the convective 
elements and has, to some extent, been borne out 
by calculations (Taylor 1973, Weir 1975). Some 
justification for limiting the angular velocity 
to some fraction of the critical velocity is also 
provided by these calculations, though the choice 
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of 0*9 for the limiting angular velocity is 

entirely arbitrary. Actually, constant specific 

■ 2 

angular momentum implies that oi “ = constant , 

where o7 is the distance from the rotation axis. 

However, within the restriction that j; be constant 

2 

on spherical surfaces j^r = constant is the closest 
possible approximation. 

V . THE EVOLUTION OF A 7M^ STAR : RESULTS 

The Case 1 sequence was terminated during 

helium burning in the core. At this point it 
was clear that the effects of rotation on the 
post -ViS evolution are negligible if solid body 
rotation throughout ti'3 star is maintained during 
these stages. Unless coupling mechanisms , such 
as magnetic fields, which do not directly depend 
on circulation currents are stronger than one 
would normally expect , solid body rotation throughout 
the star will probably not be maintained. In any 
case, the numerous previous calculations for non- 
rotating stars adequately describe the Case 1 
evolution so Case 1 will not be discussed further. 
Cases 2-4 were continued until critical velocities 
were encountered at some point in the star. In 



all three cases, this occurred after helium exhaustion 
in the core and prior to carbon ignition. 

Time Scales 

One of the effects of rotation is to lengthen 
the time scales for the evolution. Table 2 gives 
the ages of the models at a number of stages for 
Case 0 and Cases 2-4. Hydrogen and helium exhaustion 
are defined to occur when the abundances of the 
respective nuclei drop below the minimum abundance 

_4 

(10 ) allowed for in the code. Helium ignition 

is defined to occur at the first appearance of 
a convective core due to helium burning. Some 
helium burning , as evidenced by a slight decrease 
in the helium abundance, generally occurs before 
this but the first appearance of a convective 
core provides a convenient point at which to 
compare the models in the various sequences. 

Subsequent to helium exhaustion, the hydrogen 
burning shell is extinguished and the convective 
envelope intrudes into the hydrogen exhausted 
core, reducing the size of the. core. This 
stage and the stage at wt. ich the maximum depth 
of the convective envelope is reached and'the 
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hydrogen burning shell reignited were chosen as 
two more points at which to compare the sequences. 
Below the ages at each stage, the tf-ines required 
to evolve from the previous stage are given. 
Finally, for Cases 2-4, this time is compared to 
the time required by the non-rotating (Case 0) 
models. In Case 4, critical velocities were 
reached prior to the envelope intrusion stage. 
Inspection of Table 2 shows that the time scales 
for the evolution are not grossly affected by 
rotation. E.xcept for Case 4 , the increases in 
the time scales prior to helium exhaustion are 
less tlian . The 20% increase in tlie duration 
of heliu.m burning in Case 4 is larger than the 
present uncertainty in Che time scales for non- 
rotating models but it is also strongly dependent 
on the limiting angular velocity allowed in the 
convective core. 

Loops ivc the U-R Diagram 

Figure 2 shows tnc trac'.;s of the seque.aces 
for Cases 2-4 in the H-R diagram. The points 
at which the stages listed in Table 2 occur are 

indicated by the letter symbols given in Table 2. 
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For comparison, the Case 0 track is shown with 
dashed lines, Here again, the effects of 
rotation are strongest for Case 4. There appears 
to be three loops for this case. The first 
(shortest) loop is probably due to a problem in 
locating the outer boundary of the convective core. 

An underestimate of the size of the convective core 
leads to an overestimate of the decrease in the 
helium abundance in the core. This decreases the 
helium burning energy generation rate and causes 
the sequence to swing to the right in the H-R 
diagram, prematurely ending the loop. In subsequent 
models, the core grows larger again, the helium 
abundance in the core increases, and a larger loop 
follows, The third (highest luminosity) loop 
occurs during the early helium shell burning phase. 
The evolution through these stages (for Case 4 
only) was later recalculated using shorter time - - 
steps and a finer criterion for the size of the 
convective core. In this calculation, the first 
and second loops merged into a single loop i-esembling 
the loop for the Case 0 sequence , though extending 
to slightly higher effective temperatures. However, 
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the loop which occurs during the helium shell 
burning phase still appeared, indicating that it 
is real and that multiple loops can result from 

rotation. . . , 

The results of OfT and Meyer-Hofmeister 
also indicate that rotation can increase the size 
and number of loops daring these stages . One 
of the differences between our sequences and those 
of KMT and Meyer-Hofmeister is that their models 
begin v/ith critical velocities on the ZAMS whereas 
our ZAMS models are rotating about one -half as 
fast (at the surface). With respect to redistribution 
of the angular momentum, our Case 2 is equivalent 
to the Case b of KMT. In Case ac of KMT and the 
models of Meyer-Hofmeister , local conservation of 
angular momentum is assumed in all radiative 
regions and convective regions are treated as in 
Case 2. Meyer-Hofmeister has emphasized that the 
effects of rotation may be impoi’tant for the 
statistics of cepheids, specifically by lowering 
the smallest mass for which loops crossing the 
instability strip occur. This would presumably 
shift the peak in the number vs. mass distribution 
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of cepheids to a lower* mass. Meyer-Hofmeister 
estimates that the peak will.be shifted by about 
one solar mass below the peak for non -rotating 
stars. There are , however, several points which 
should be made in this regard. First, the loops 
in the H-R diagram represent one of the more 
uncertain features of stellar evolution theory, 
even for non-rotating stars. As is discussed 
by Meyer-Hofmeister, the number and length of the 
loops are very sensitive to the handling of the 
convective core boundary during hydrogen burning. 
Thus, the effects (discussed in Section II) of 
rotation on the criterion for convection should 
be taken into account. This was not done at all 
in the calculations of KMT and Meyer-Hofmeister 
and only approximately (one of the effects) 
in the present investigation. Second, our models, 
which represent typical stars (rather than rapid 
rotators), do not show significant effects of 
rotation in the H-R diagram, except in Case 4. This 
is admittedly an extreme (though possible) case. 
Finally, it is likely that stars w’hich begin their 
evolution with near-critical velocities will lose 
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a considerable fraction of the angular momentum in 
their envelopes prior to helium burning. Foi- 
example , equatorial mass shedding during the Be 
phase has been shown by Strit,tmattei', Robertson, 
and Faulkner (1970) to be a very efficient means 
of removing angular momentum. Kraft (1966) 
has shown that the observed rotational velocities 
of cepheids are consistent with the loss of a 
considerable fraction of the MS angular momentum. 

It appears , therefore , that furtner investigations 
are requii-ed before any sti'ong conclusions can 
be reached about the effects of rotation on t.he 
statistics of cepheids. 

One fui'ther comment should be made about 
Figure 2 . „ At point E, the giant branches of our 
sequences show a pi'onounced bending to the right. 
This is the point at which intrusion of the con- 
vective envelope into the hydrogen exhausted core 
begins. This mixes material from the core into 
the envelope and lowers the hydrogen abundance 
there. To take this into account in the envelopes, 
the opacities from the tables (which applied to 
the ZAMS composition) were reduced as if the 
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opacity was due to pure electron scattering in an 

ionized medium. Since much of the mass in the 

envelopes during this stage lies above the hydrogen 

ionization zone, the affect of the reduction of the 

hydrogen abundance in the envelope was overestimated 

by this method. 

Conditions at the Center 

Figure 3 shows the paths of the sequences 

for Cases 2-4 in the (log central density, log 

central temperature )“plane. Again, the stages 

referred to in Table 2 are indicated by capital 

letters and the Case 0 sequence is shown by 

dashed lines. Figures 4 through 6 show the 

angular velocities (solid lines) of the inner most 

mass zone as functions of the model numbers. The 

dashed lines show a)A> for the same zone. The 

cr 

letters at the tops of the figures refer to the 
evolutionai'y stages. The discontinuous nature of 
the angular velocity curves for Cases 2 and 4 are 
due to the finite zoning of the models. For 
instance, after helium exhaustion (point D) , the 
helium exhausted core grows by finite increments 
as the helium burning shell moves outward. Whenever 
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the helium abundance in a shell reaches zero and 
; that shell becomes part of the core (which is 
rotating as a solid body), the shell’s angular 
momentum and moment of inertia are included 
in computing the angular velocity of the core in 
the next model. This causes a sudden drop in the 
core angular velocity. 

From Figure 3, it is apparent that, for Cases 
2 and 3, rotation does nn+ significantly affect 
the central conditions until after helium exhaustion 
(point D). For Case 4, however, there are 
noticeable effects during both the hydrogen and 
helium burning stages. Figure 6 shows that the 
angular velocity of the inner shell is at the 
limiting velocity during these stages and, 
therefore , the magnitudes of these effects depend 
on the chosen limiting velocity. After convection 
in the core dies out and before the fuel is completely 
exhausted, there is a shox’t period when the core 
rotates with local conservation of angulai' momentum 
and with no limiting velocity. After the hydrogen 
burning stage , the angular velocity reaches 99% 
of , during this period. 
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After helium burning, rotation becomes more 

and more important in all cases as the core contracts 

and the angular velocity increases rapidly. In 

general, rotation causes the temperature , at a 

given density, to be lower in the rotating models 

than in the non-rotating models. Thus, carbon 

ignition will become very difficult to achieve 

unless some mechanism couples the core to the envelope 

In view of the large gradient in the mean molecular 

weight, across the hydrogen burning shell, it is 

difficult to see how this coupling will occur. 

(However, KOT have suggested some possible mechanisms. 

The evolution of the angular velocity distribution 

After hydi'ogen burning, the contractions 

and expansions of various parts of the star plus 

the assumed rotation laws produce strong differential 

rotation. This is illustrated in figures 7 to 9, 

which show the angular velocities and y. /■\^ as . 

'cr 

functions of the mass coordinate foi' Cases 2-4, 
respectively. F-ach figure shows models at stages 
A, C, and D and a inodel close to the end of the 
sequence. A general feature of all the cases 
is the spin-up of the core. Critical velocities 


were reached in each case soon alter the last 
model shown in the figures. The critical velocities 
were reached, not at the center, but rather, at 
the outer edge of the hydx'ogen exhausted (Case 3) 

Or helium exhausted (Cases 2 and 4) core. This 
situation can be seen, for Cases 2 and 4 , as a 
direct result of imposing solid body rotation in 
the chemically homogeneous core. For Case 3, the 
occurrence of critical velocities at the -outer 
edge of the core can be understood by considering 
the structure of a nearly isothermal, degenerate 
core. For a non-re lativistic , degenerate core, 
the density distribution is very similar to that 
of a n = 1.5 polytrope (Chandrasekhar 1939). The 
density distribution in the convective core of 
the ZAJIS model can also be approximated by a n = 1. 5 
polytrope. For local conservation of anguJar 
momentum, the evolution of the angular velocity 
distribution is determined by the evolution of 
the density distribution. Thus the angular 
velocity distribution in the degenerate core of 
the final model is very similai* to that of the 
ZAMS model, i.e,. flat, though the value of the 
angular velocity in the final model is almost 
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four orders of magnitude larger. This nearly flat 

distribution causes m/iij to increase outwards 

cr 

from the center. 


In all of our rotation laws, we have ignored 
the possibility that instabilities can overcome 
gradients in the mean molecular weight (a -barriers) 
F.owever, for large enough gradients in the angular 
velocity and small enough gradients in the mean 
molecular weight, instabilities can occur (cf. 

Zahn 1974). In the region of transition from 
the carbon-oxygen core to the helium-rich zone , 
the gradient in the mean molecular weight is 
fairly small and, in Cases 2 and 4, the gradients 
in the angular velocity are very large. If we 
had included such instabilities, redistribution of 
the chemical composition, as well as the angular 
momentum, would have occurred. The effect this 
would have on the models is difficult to predict 
without detailed calculations since the effect 
would depend on the interaction of several factors: 
the degree of mixing produced by the instability, 
the effect of the mixing on the helium burning 
shell, the timescale for re-establishing something 
like solid body rotation in the core, etc. 
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VI . DISCUSSION AND CONCLUSIONS 

The formulation we have used for the stellar 
structure equations assumes that the total potential 
is conservative. Our rotation laws (constant n on 
equipotential surfaces) are rot conservative, 
except in regions with solid-body rotation. The 
result cf such non-conservative i-otation laws 
will be that p and T are not strictly constant on 
equipotential surfaces. Although vve have not 
explored in detail the effects this would have 
on the models, we may hope that any such effects 
will be smoothed by circulation currents which 
generally operate very efficiently along equipotcntials . 
The numerical tests presented in Section III 
indicate that cur models agree well with models 
employing conservative potentials. 

From the results presented in Section V, 
we may draw several conclusions about the role 
of rotation in the evolution of a star leaving 
the MS with an average rotation rate: 

1. rotation does not severely a.^fect the total 
time scales for evolution because rotation 
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does not become important until the later 
stages , where the time scales are intrinsically 
very short ; 

2. The effects of rotation on the loops in the 
H-R diagram are still very uncertain and 
depend critically on the initial rotation rate 
and the distribution of the angular momentum 
with the star; and 

3. unless some mechanism is capable of transporting 
angular momentum across u-barriers, critical 
velocities are reached prior to carbon 
detonation, regardless of the rotation law 
within chemically homogeneous regions. 

It appears that secular instabilities of 
the Goldreicii -Schubert type may be able to overcome 
the weak u-barrier at the outer edge of the carbon- 
oxygen core and redxstr xbute the angular mornGntum.' 
Since such instabilities will cause mixing of 
helium into the carbon-oxygen core, further model 
calculations are necessary to determine what 
effects this will have on the evolution. Such 
calculations will be the subject of the next paper 
ir. this series. If it turns out that it is 
necessary to transport angular momentum across the 
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strong a -barrier in the region of the hydrogen 
burning shell , more extreme redistribution 
mechanisms , such as those suggested by KNJT, may 
have to be explored. Of coure , the role of 
magnetic fields in these stages is still largely 
an open question. 

Although the results presented in this paper 
are not complete enough to determine exactly 

what the ' effects of rotation are on the post - 

helium burning stages, it is clear that rotation 
will play a dominant lolc. The evolution of a 
rotating star will be significantly different 
from tliat indicated by models of non-rotating 
stars. Aside from the significant effect of the 
reduction of the effective gravity , rotationally- 
induced mi.xing in regions which are stable to 
purely thermal convection will change the size 
arid chemical composition of the core and the 
behavior of t-he sTriell sources. This may tr.tally 
alter the evolution in the subsequent stages. 
Finally, we emphasize that the above considerations 
apply to stars which leave the main sequence 
with moderate rotational velocities, as well 
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as stars which may be considered "rapid -rota tors" 
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where r is the radius oi the equiootentia 1 surface 
at the angle defined such that P^Ccos = 0, 

P,, is the second-order Legendre polynomial , and 
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Y. is the axisymmetrlc tesseral harmonic relating 
3 

r to r^ on a given equlpotentlal surface: 

r(r^, 0) = r^ [l + ^ 0)j. (A4) 

■ ■ j. ■ ■ . 

Consistent with the approximation for mentioned 
above, Y. is given by , 



2 , 2 



(cos 6 ) 


(A5a) 


and 


Y. = O for j jt 2, 

s) ■ 


(A5b) 


where Tig is the logarithmic derivative of Yg with 
respect to r^ : 




(A6) 


This quantity can be evaluated by integration of 
Radau's equation (cf. Kopal 1959) - 

+ ’Ij ‘’Ij - » - 

j(j +1) (A7) 

with j = 2 and the boundary condition 71^(0) = j - 2. 
Using equation (A6) we can evaluate the derivative 
in equation (A3) as follows : 
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lr-(r^^Yo)= 5 




A ^'^9 



or, 

^ (r^^ Yg) = r/ (5 + ^ . (A8) 

At this point, we can write the total potential 
as 





(cos 9 )q 




2 

■ju 


^ + ‘!2 
2 + ilg 



+ ^ ^ ju^ Pg (cos 9) + 

const . , (A9) 






where all the quantities within the integral are 
to be evaluated at the point referred to by the 
dummy -variable-of -integration . 


If we define 



2 ^ 2 
(u r^ 



(AlO) 


then the equation of an equipotential surface Is , 
from equations (A4) and (A5) , 


r = [ 1 - A ?2 (cos 9)]. (All) 

In order to relate r„ to r. , we evaluate the 

o , ijf ’ 

volume integral from r = o to r given by equation 
(A9) . This gives 

% = ^ V (A12) 

so, from equation (1) 

•3 0 . 0 1 ^/3 

— o t l - • (A13) 

In practice, one knows r and wants to find r 

■f o 

Since A is a function of r^ , this must be done by 
an iterative procedure . - 

Since the local effective gravity is given by 



r'i M'i 

\r aey ’ 


(A14) 
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g can be found by differentiation of equation 
(A7) . The integral in equation (AB) and its 
derivatives must be evaluated numerically. Once 
g and g~^ are known for a set of points on an 
equipotentia 1 surface , S g and S, g”^ can be 
found from equation (2) by numerically integrating 
over (5 . ■ . " . . • 

In the procedure we have used, V 12 is not 
evaluated by direct integration of equation (A7). 
Instead, the ratio is evaluated at the 

point p and is assigned the value appropriate 
for the surface of a polytrope with the same value 
of The values of (as a function of the 

poly tropic index) wei’e taken from Kopal (1959) 
and the values of p/p^ from Chandrasekhar (1939) , 
Test models computed with direct integration of 
eqviation (A 7) showed that this approximation has 
no discernible effect on the models. At this 
point, all the quantities needed in equations (10) 
and (16) to evaluate fp and f.^, are known. 
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APPENDIX B 


IMPLEMENTATION IN THE STELLAR 
EVOLUTION CODE 

ENVELOPES 

The Paczynski code uses a separate program 
to calculate a grid of envelopes in the (log L, 
log Tq) - plane , where T^ is the temperature of 
the outermost shell in the envelope (for a 
complete description of this program, see 
Paczynski 1969) . These envelopes are then used 
to provide outer boundary conditions for the 
interior (Henyey) program. Because the envelopes 
are constructed before the interior structure is 
known, the integral in equation (A9) cannot be 
directly evaluated and the following approximation 
is used. 

Let B denote the integral in equation (A9). 

We evaluate B, for the envelope only, as if the 
interior is a polytrope rotating as a solid body. 
Then uj and T |2 are constant and can be taken outside 
of the integral. For a poly trope of index n, 
p/o^ = S^> where q is the Lane-Emden function 
(not to be confused with the angle 9 ).. If we 
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introduce the variable § = (tq/R) where R is 
the radius of the polytrope and is the first 
zero of the Lane-Emden function, then 

„ VR A7 ® * "2 r? e” 5'® j,/ % 

■ Jo ® V ‘ ’ 

The interior mass in a polytrope obeys the equation 


•5 (- 5^ !!)■ 


® 5-+-Ti" j t-'^/Fe 

X iS O 


d§'. (S3) 


Equating R with r^ and rewriting the above equation 


^ "^2 r^i 


IF 7T 


.h . ,4 

j.vK - n - 


yields an equation which contains only quantities 
related to polytropes on the right-hand-side. 

The solutions of the Lane-Eraden equation given 
by Comrie (1932) were used to evaluate equation 
(B4) for n =1.5, 2, 2.5,.,., 5 and the results 
put into the envelope program as a table. The 
same procedure for interpolation in as used 
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for finding iTg (see Appendix A) is then used to 
find 

This procedure requires that d and be 
known beforehand, it is assumed that ju is constant 
throughout the envelope and the 2-dimensional 
grid of envelopes if replaced by a 4-dimensional 
grid in (log L, log T^, ji , p^). In practice, 
the envelope program is used as a subroutine 
of the Henyey program which calls for new 
envelopes whenever the models run out of the old , 
grid. In this way, the grids can be kept small 
without having to frequently stop the run to 
calculate new envelopes. 

INTERIORS 

For the interior (Henyey) part of the 

calculation, the equations in Appendix A are 

used, with no further approximations. . A flow 

chart for the program is given in Figure Bi . " ~ 

Two tests are employed for convergence of the 

rotation factors f^ and f,j,; if eicher test is 

passed, the model is assumed to be converged. 

Generally, TEST 2 is the determining factor. This 

means that f and f,., have converged to the point 
P ; * 

where further changes will not affect the model 



beyond the tolerances set by the Henyey program, 
i.e:, the model is self-consistent. TEST 1 
was put in primarily to allow the program to be 
used for non-rotating models without requiring 
two passes through the Henyey program for each 
model. For rotating models, two or three passes 
through the rotation loop are usually sufficient 
though as many as ten iterations may be required 
when the radius of the star is rapidly changing. 
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♦Abbreviations : 

SB - solid body rotation and overall conservation of 
angular momentum; 

LC - local conservation of angular momentum; and 
2 

CA - a'l* =• constant (but r i 0.9 and overall 

conservation of angular momentum (see text). 
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TABLE 2 

EVOLUTION TIME SCALES* 


STAGE 

SYMBOL 

CASE 0 

CASE 2 

CASE 3 

CASE 4 

ZAMS 

A 

0.0 

0.0 

0.0 

0,0 

H-exhaustion 

B 

t(A-B) 

32.861 

32.861 

33.208 

33.208 

1.011 

33.676 

33.676 
1.025 

33.706 

33.706 

1.026 

He-ignition 

C 

-r (B-C) 
'/'o 

35.392 

2.531 

35.803 

2.594 

1.025 

36.311 

2.635 

1.041, 

36.475 

2.769 

1.094 

He -exhaustion 

D 

t(C-D) 

^/■o 

50.806 

15.415 

51.691 

15,888 

1.031 

52.233 

15.922 

1.033 

55.049 

18.573 

1.205 

Envelope 

intrusion 

E 

51 .767 
0,961 

52.842 

1.151 

1.198 

53.506 

1.273 

1.325 

— 

H-reignition 

F 

T (E-F) 
tAq 

51.870 

0.103 

52.962 

0.120 

1.165 

53.645 

1,139 

1.35Q 



6 

♦Time scales are in units of 10 years . 
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FIG. 1 


FIG. 2 


FIG. 3 


FIGURE CAPTIONS 

- The effects of differential rotation on a 
30M^ ZA&!S model as predicted by the KT 
method (dashed line) and by the SCF 
calculations of Bodcnheimei* (1971, symbols) . 

The zero subscript refers to the non- 
rotating model, which is compared to 
rotating models with different values of 
the total angulax momentum (J) . For the KT 
models , the angular velocity distribution was 
spherical with constant (in the equatorial 

plane) throughout the star. Values of 

are given at the top. For the SCF models, 
the different symbols refer to different 
angular velocity distributions and are the 
same as used in Figure 5 of Bodenheimer (1971) , 

- The H-R diagrams for Cases 2-4. The stages 
referred to in Table 2 are indicated by 
captial letters. The Case 0 track is 
indicated by dashed lines. 

- The paths in the (log log T^)-plane 

for Cases 2-4. Stages are indicated by 
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FIG. 


FIG. 

FIG. 

FIG. 


FIG. 

FIG, 

FIG. 


capital letters. The trajectory of the 
Case 0 sequence is indicated by dashed 
lines . 

4 - The angular velocity at the inner most 

mass zone as a function of MODEL NUMBER 
for Case 2. The angular velocity in the 
equatorial plane is shown by a solid line 
and o)/?' by a dashed line. The numbers 
are rounded to two decimal places. 

5 - Same as Figure 4 , but for Case 3 

6 - Same as Figure 4, but for Case 4 

7 - The angular velocity distribution through- 

out the models at four stages in the Case 
2 sequence. The angular velocities are 

indicated by solid lines and by 

CX*' 

dashed lines , 

8 - Same as Figure 7, but for Case 3 

9 - Same as Figure 8, but for Case 4 

B1 - Flow chart for the rotating model program. 
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